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'. ABSTRACT 



1 INTRODUCTION 

The intergalactic medium (IGM) contains most of the 
baryons in the Universe and it provides the fuel for galaxies 
to form stars in which metals are produced. In turn, super- 
novae and galactic winds enrich the IGM with metals, while 
stars and active galactic nuclei (AGN) emit UV photons. 
This interplay between the IGM and galaxies, mediated by 
metal cooling in the presence of UV, regulates the forma- 
tion of stars in the universe. The evolution and enrichment 
history of the IGM provides a record of this interplay. 

Observations of metal absorption lines (e.g., C III , C 
IV, Si III, S IV and O VI) in quasar spectra show that 



A study of metal enrichment of the intergalactic medium (IGM) using a series of 
smooth particle hydrodynamics (SPH) simulations is presented, employing models for 
metal cooling and the turbulent diffusion of metals and thermal energy. An adiabatic 
feedback mechanism was adopted where gas cooling was prevented on the timescale of 
supernova bubble expansion to generate galactic winds without explicit wind particles. 
The simulations produced a cosmic star formation history (SFH) that is broadly con- 
sistent with observations until z ~ 0.5, and a steady evolution of the universal neutral 
hydrogen fraction (rini) that compares reasonably well with observations. The evolu- 
tion of the mass and metallicities in stars and various gas phases was investigated. At 
z=0, about 40% of the baryons are in the warm-hot intergalactic medium (WHIM), 
but most metals (80%-90%) are locked in stars. At higher redshifts the proportion of 
metals in the IGM is higher due to more efficient loss from galaxies. The results also 
indicate that IGM metals primarily reside in the WHIM throughout cosmic history, 
which differs from simulations with hydrodynamically decoupled explicit winds. The 
metallicity of the WHIM lies between O.OI and O.I solar with a slight decrease at lower 
redshifts. The metallicity evolution of the gas inside galaxies are broadly consistent 
with observations, but the diffuse IGM is under enriched at z ^ 2.5. 

Galactic winds most efficiently enrich the IGM for halos in the intermediate mass 
range IO^^Mq - 10^^ Mq. At the low mass end gas is prevented from accreting onto 
halos and has very low metallicities. At the high mass end, the fraction of halo baryons 
escaped as winds declines along with the decline of stellar mass fraction of the galaxies. 
This is likely because of the decrease in star formation activity and decrease in wind 
escape efficiency. Metals enhance cooling which allows WHIM gas to cool onto galaxies 
and increases star formation. Metal diffusion allows winds to mix prior to escape, 
decreasing the IGM metal content in favour of gas within galactic halos and star 
forming gas. Diffusion significantly increases the amount of gas with low metallicities 
and changes the density-metallicity relation. 

ISimcoe et aLlliooi ). There is evidence for enrichment ex- 
tending back to z > 5 (|Pettini et al.1 l2003l : ISimcod |2006| ) 
though metallicitie s do not evolve much from z=4 to z=2 
(Sch ave et al.|[20o3 ). Since metals are created in stars inside 
galaxies, those in the IGM must have escaped from galax- 
ies. Exactly how this occurs is unclear. It typically assumed 
that galactic winds, driv en by star formation, dominate IGM 
enrichment (lAguirre fc S chave 2002). Winds can be either 
launched as the ejecta of a large number of co- existing su- 
pernova (SN) explosions (|Heckman et al.lll990l ). driven by 
the injection of momentum by SN and stella r winds or by 
radia tion pressure from starbursts and AGN urrav et alJ 
2005j). Observation s of galactic winds (e.g., Hccknian 20oJ|; 
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Pettini et al.ll200ll ) have found that galactic wind velocities 



range from hundreds to thousands of km/s and the mass loss 
rates are comparable to star formation rates. Wind material 
has a complex, multiphase structure but most metals are 
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expec ted to be entrained in the hot phase (|Veilleux et alj 
I2OO5I . and references therein). 

Although detailed hydrodynamical simulations 
of the interstellar medium (ISM) in galaxies (e.g, , 
Mac Low fc Ferraral Il999l ; [Strickland fc Stevens! I2OO0I : 
WiUiams fc DvsonI I2OO2I ) have been able to generate 



galactic outflows and have explored various properties of 
winds and the gas dynamics in different phases, current 
cosmological simulations lack the resolution to launch or 
track winds directly. Hot, low density SN bubbles are 
unresolved in such simulations which initially led to an 
overcooling problem that p roduced unrealistically con - 
centrated simulated galaxies l|Navarro fc Steinmeta [T997h . 
As a result, various "subgrid" stellar feedback and wind 
models have emerged. These models serve two functions: 
to regulate star formation and the properties of the ISM 
and to redistribute gas (and newly formed metals) both 
within and into the environment around galaxies. There are 
three main approaches, energetic feedback, kinetic feedback 
and modifications to the effective equation of state which 
behaves similarly to an increased effective pressure. 

Energetic feedback in its simplest form involves simply 
adding the stellar feed back as thermal e nergy, but this suf- 
fers f rom overcool ing ("Katz et al."l996^. Kinetic feedback 
(ej 



Navarro fc W hite 1993; Springcl & Hcrnquist 



iQppenheimer fc Pavel 12009 : IDaUa Vecchia fc Schavd l20oi ') 
converts part of the SN energy into kinetic energy in the 
gas. The effectiveness of kinetic feedback is strongly depen- 
dent on the resolution and hydrodynamic method. 

ISpringel fc Hernguistl |2003) argued that regulated star 
formation creates an effective pressure in the ISM and this 
was modelled directly in the GADGET code as part of a 
recipe for regulated star formation. This approach leads 
to a strongly hydrodynamically-coupled, multiphase ISM 
that does not naturally produce galactic outflows. To com- 
bat this the authors added a "superwind" model where 
fluid elements in the star forming region are ejected at 
flxed speed and are als o hydrodynamically decouple d until 
they leaves the galaxy. lOppenheimer fc Davd (12006') modi- 
fled the model in a manner referred to as the momentum- 
driven wind scenario so that the velocity of the wind and 
the mass loading factor were related to the velocity dis- 
persion of the host galaxy. The GADGET code with super- 
wind feedback prescriptions has been widely used in var- 
ious problems such as damped Lyman-a (DLA) absorbers 
l|Nagamine et al ]|2004i|l and the enrichment of the IGM at 
high and low redshifts (joppenheimer fc Pavel l2006l . Hooi). 
According to these works, superwind feedback is essential 
to suppress overproduction of stars in galaxies and to repro- 
duce the cosmic SFH at high redshift. It also increases the lo- 
cal fraction of the warm-hot intergalactic medium (WHIM) 
to a sufficient percentage (40% to 50%) to account for the 
"missing baryons" at z = 0. Although aspects of the model 
compare well with observations, some components do not. 
For example, the feedback may eject a large amount of cool 
gas from the galactic disks, which r esults in a low neutral hy - 
drogen mass density Hhi at z < 2 ijNagamine et al.ll2004al ). 
Also, the interaction between wind s and the ISM is usually 
not m odeled in these simulations. iDalla Vecchia fc Schavd 
l|2008l ) found that the ISM plays an important role in regu- 
lating the amount of wind that escapes and the morphology 
of the galaxies. In their model winds are not hydrodynam- 



ically decoupled, which naturally allows for variable mass 
loading. 

A refined version of energetic feedback is adiabatic 
feedback which treats the overcooling problem by inhibit- 
ing gas from cooling until the hot SN b u bbles can be re- 
solved (e.g., Thacke r fc Couchman* '2000*; ^Kav et ah I I2OO2I : 
ISommer-Larsen et a l. 2003; Stinson ct al. 2006). This is the 
approach used in this work. The pressure of the hot gas ac- 
celerates the ISM to generate winds. In the high resolution 
limit, this method approaches direct ISM modeling. Though 
energetic feedback is often referred to as supernova feedback, 
it can be used to model several types of stellar feedback such 
as winds and locally deposited radiation energy. The essen- 
tial quantity is the energy injection rate as a function of the 
mass in stars and the current age of the stellar population. 

iTheuns et al] (|2002l ') used the iKav et al.1 l|2002l ') adia- 
batic feedback model that turned off cooling for 10 Myrs 
for the feedback gas in their cosmological simulations. They 
found enough metals were carried by strong winds to pro- 
duce C IV ab s orptio n lines that agreed with observations. 
lAguirre et al.1 l|2005h used the same simulation to compare 
the optical depth of C IV and C HI absorption lines from 
simulations with observations. 

The properties of the enriched IGM are not only affected 
by winds but also by gas cooling. Winds can enrich galac- 
tic halos and the IGM so that metal cooling significantly 
increases the cooling rates. lAguirre et al.l |200g) found that 
their simulated metal enriched gas was too hot (10^ ~ 10^ 
K) and suggested that a lack of metal cooling was re- 
sponsible for discr epancies between simu l ated a nd observed 
C IV absorption. lOppenheimer fc Pavel ||2006D found bet- 
ter agreement when they included the Sutherland fc Popital 
(199 j) metal cooling model. The same model was used 
by Choi fc Nagaminef (|2009| ) , who investigated the effect of 
metal cooling on galaxy growth and found that it increases 
the local star formi ng efficiency and enha n ces ac cretion onto 
galaxies. However, ISutherland fc Popital l| 19931 ) did not in- 
clude photoionization due to a ultraviolet (UV) radiation 
background which strongly affects the ionization states of 
metal species and changes the cooling rates. This is investi- 
gated in detail in the current work and also in |w iersma et al.l 
|2009i). 

Another important aspect of metal enrichment is the 
mixing of metals between the wind and the surrounding 
gas. The interstellar medium is highly turbulent and SN 
explos ions are likely to be a m ajor driver of the turbu- 
lence (|Mac Low fc Klessen|[2003 '). In addition, large velocity 
shear (such as between a wind and a gaseous halo) naturally 
generates turbulence and mixing. Turbulent mixing redis- 
tributes metals and thermal energy between the wind fluid 
and the ambient gas. This changes the metallicity, temper- 
ature and future evolution of the gas. While metal mixing 
is expected in strong outflows, it is still unclea r how mixing 
impact the IGM. For example, observations bv lSchave et al.l 
(|2007l ) found compact (~ 100 pc), transient C IV absorbers 
that are highly enriched, suggesting poor chemical mixing at 
small scales. These absorbers were interpreted as enriched 
clumpy medium embedded within hot galactic wind fluids. If 
velocity shear is the major mechanism for turbulent mixing 
between winds and the surroundings, then this poor mixing 
could be explained if the clouds are carried by hot winds 
at the same speed. However to investigate this in detail. 
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one must resolve wind structures, which is beyond current 
cosmological simulations. In this work we will focus on sub- 
grid turbulent mixing models in cosmological context. In 
SPH simulations (which represent the majority of work in 
this area), the fluid is modeled by discrete particles. This 
implies that newly injected metals are locked into specific 
particles. For example, it was found that the distribution 
of metals from SPH simu lations is too inhom ogeneous com- 
pared with observations (lAguirre et al.ll2005l). To asse ss the 
potential importance of mixing, IWiersma etall (|2009bh used 
SPH-smoothed metallicities and compared it with conven- 
tional particle metallicities, and found smoothing is able to 
generate significantly more material with low metallicities. 
This approach cannot capture the spread of metals over time 
with its impact on cooling and the thermal history of the gas. 
Directly modeling the turbulent ISM within a cosmological 
simulation is far beyond curr ent ca pabilities. We employed 
a variant of the fSmagorinsk^ l| 19631 ) subgrid turbulent diffu- 
sion model, in which unresolved turbulent mixing is treated 
as a shear-dependent diffusion term. Metal cooling was cal- 
culated based on the diffused metals so that its non-local 
effects could be investigated. 

In this work, we present an analysis of a series of SPH 
cosmological simulations that incorporated adiabatic stellar 
feedback, detailed metal cooling and turbulent mixing to 
study the evolution and enrichment of the IGM. The feed- 
back model was kept si mple, following the a diabatic stel- 
lar feedback approach of IStinson et al T (|2006D . This model 
has b een calibrated via num erous galaxy formation studies 
(e.g., iGovernato et al.]|2007t ). No additional wind prescrip- 
tions were used. With this approach, outflows arise from 
stellar feedback within the ISM and there is no distinction 
between the feedback that regulates star formation and that 
which drives galactic outflows. Thus this work establishes 
a baseline for the effectiveness of moderate stellar feedback 
coupled with key physical process absent from other work to 
reproduce the properties of the IGM. These results may be 
compared with explicit wind models. A further goal of this 
paper is to separately quantify the impact of metal cooling 
and turbulent mixing on the SFH, the global properties and 
the evolution of the IGM and its enrichment. In this first 
paper of a series, we present general results. The properties 
of specific metal absorbers in simulated quasar spectra will 
be presented in a second paper. 

This paper is organized as follows. Section [5] describes 
the models for cosmological hydrodynamics, star forma- 
tion, supernova feedback, metal cooling and metal diffu- 
sion. Section [3] examines the cosmic SFH, global H I frac- 
tion and Ly-a decrement in order to calibrate our models. 
Section |4] focuses on the evolution of the baryonic mass, 
metal fractions and metallicities in stars and different gas 
phases. We compare those results to the observed metal 
fractions and metallicities at different epochs, and with the 
simulations using diff e rent s ubgrid feedback model s from 
lOppenheimer fc Pavel (l2006l ). iDave fc Qppenheimerl (|2007l ) 
and Wiersma et al.l I 2009bl ). In section[5]we analyze the dis- 
tribution of mass and metallicity in the density-temperature 
phase diagram at z = 0. In section[6]we characterize our wind 
efflciency as a function of galaxy mass to obtain a better un- 
derstanding of how different phases of the IGM get enriched. 
Where relevant, we have included detailed analysis of the ef- 
fects of metal cooling and diffusion, and comparisons with 



Table 1. List of Simulations 



Name 


Size (Mpc) 


Np 


Metal Cooling 


Diffusion 


mcd_40_256 


40.0 


2 X 256^ 


Yes 


Yes 


nmc_40.256 


40.0 


2 X 256^ 


No 


Yes 


nmd_40.256 


40.0 


2 X 256=* 


Yes 


No 


mcd_45.256 


45.6 


2 X 256=* 


Yes 


Yes 


mcd_45_512 


45.6 


2 X 5123 


Yes 


Yes 



observations. In the final section [T] we summarize and dis- 
cuss the broader implications. 



2 METHOD AND SIMULATIONS 

Table [T] lists the simulations used in this study. The cosmo- 
logical parameters are {Qrn.,^A,^b, h, as, n) — (0.279, 0.721, 
0.0462, 0.701, 0.769, 1.0 ), consistent with the WMAP 5-year 
results (jKomatsu et al The first three simulations in 

the table share the same initial condition but have different 
subgrid models for metal cooling and diffusion. The parti- 
cle mass for gas is rrig — 2.4 x IO^Mq and for dark matter 
is rud = 1.2 X 10** Mq. When discussing he effects of metal 
cooling and diffusion, we will usually refer to the first simu- 
lation (mcd_40_256) "the reference run". The last two sim- 
ulations were designed for a convergence study. They have 
the same initial condition but differ from the first three. The 
high resolution case has particle masses md = 3.6 x IO'^Mq 
and TUg — 4.5 x lO'' Mq. All the simulations have the same 
star formation, metal production and SN feedback models 
described below. The simulations were run to redshift zero, 
except for the "mcd_45_512"case which was stopped at z = 
2. 

The simula tions were evolved u sing the parallel SPH 
code GASOLINE IWadslev et aLlliooil ). gasoline solves the 
equations of hydrodynamics and includes radiative cooling. 
Gravity is calculated for each particle using a binary tree 
elements that span at most 8 — 0.7 of the size of the tree 
element's distance from the particle. 

The star formation and feedback recip es used were 
the "b lastwave model" described in detail in IStinson et al.l 
(l2006h . and they may be summarized as follows. Gas par- 



ticles must be dense {n„ 



0.1 cm'-') and cool (Tn, 



15,000 K) to form stars. A subset of the particles that pass 
these criteria are stochastically selected to form stars based 
on the commonly used star formation equation. 



Mgas 



(1) 



dt tdy n 

where is mass of stars created, c* is a constant star 
formation efficiency factor, Mgaa is the mass of gas creating 
the star, dt is how often star formation is calculated (1 Myr 
in all of the simulations described in this paper) and tdyn 
is the gas dynamical time. The constant parameter, c* was 
set to 0.05 s o that a simulate d isolated mod el Milky Way 
matches the iKennicuttI (|l998l ) Schmidt Law l|Stinson et al.l 
|2006| ). Gas particles passing all the criteria form one star. 
The mass of star formed remains constant throughout the 
simulation and is set as one-third of the original gas mass 
(8.0 X 10'^ Mq for the 40 Mpc, 256^ runs and 1.5 x 10*^ Mq 
for the 512^ run). 
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At the resolution of these simulations, each star particle 
represents a large group of stars. Thus, each particle repre- 
sents a stellar popul ation covering the entire initial mass 
function presented in lKroupa et alj Star masses are 

conve rted to stellar lifetimes as described in iRaiteri et alJ 
l|l996t ). We implemented feedback from SN II, SN la and 
stellar winds. Met al enrichment from SN II and SN la fol- 
lows the model of iRaiteri et al.l (|l996l ). but metal produc- 
tion of AGB stars is not included. Stars larger than 8 Mq 
explode as SN II during the timestep that overlaps their 
stellar lifetime after their birth time. The explosion of these 
stars is t reated using the analytic m odel for blastwaves pre- 
sented in McKee fc Ostrikeil (| 19771 ) as described in detail in 
IStinson et al. I l|2006l ). Gas cooling is suppressed during the 
expansion of the SN bubble. While the blast radius and the 
cooling shutoff time are calculated using the full energy out- 
put of the supernova, less than half of that energy is trans- 
ferred to the surrounding ISM, Esn = 4 x lO^'' ergs. The 
rest of the supernova energy is assumed to be radiated away. 
The energ y is transferred to the IS M by volume weighting 
(as used in lMashchenko et al.l l|2008l l). Each affected gas par- 
ticle, with mass rrii and density pi receives a fraction of the 
SN energy and metals proportional to rmWij / pi, where Wij 
is the SPH smoothing ker n el. A s imilar method was also 
adopted in IWiersma et"ai] l|2009bl ) for the distribution of 
metals and the energy from SN la. The Blastwave model, as 
currently implemented, has a bias for more cooling suppres- 
sion at earlier epochs relative to other approaches. Numer- 
ical feedback recipes differ substantially from code to code 
and their results vary with resolution. This provided moti- 
vation for the characterization of star formation regulation 
and wind generation in these simulations that is presented 
in section [6l 

For SN II, metals produced in stars are released as the 
main sequence progenitors die and distributed to the same 
gas within the blast radius as is the supernova energy ejected 
from SN II. Iron and Oxygen are produced in SN II accord- 
ing to the a nalytic fits used in Raiteri et al.l (|l996l ) using the 
yields from IWooslev fc WeaveJ (|l995f h 



Mpe = 2.802 X IQ-^ mI'^^* 
Mo = 4.586 X lO'^Af^^-^^^ 



(2) 
(3) 



i99e) 



Feedback from SN la also follows the IRaiteri et al 
model, as described in detail in IStinson et al 
20061 ) . Radiative cooling was not disabled for SN la. 
Each SN la produ ces 0.63 Mq Iron and 0.13 Mq Oxygen 
l|Thielemann et al. 1986) and the metals are ejected into the 
nearest gas particle fo r SN la. Stellar wind feedback was 
implemented based on lKennicutt et al.l (|l994l ). and the re- 
turned mass fraction was determined using a function de- 
rived bv lWeidemamil (|l987l ). The returned gas has the same 
metallicity as the star particle. 



2.1 Metal Cooling under the UV Radiation 

CoUisionally excited metal ion species can significantly en- 
hanc e the cooling of the IGM by their line transitions (see 
e.g., [Sutherland fc Dopital l|l993l )). Metal cooling thus af- 
fects the gas temperature, the ionization states of the ob- 
servable metal species and ultimately the dynamics. It was 



not feasible to calculate metal cooling rates during the sim- 
ulation. Instead, tabulated c ooling rates were interpolate d 
during the simulations as in JOppenheimer fc Dave (2006"), 
IWiersma et~aLl (l2009bl ) and lchoi fc Nagamine (2009). Metal 
cooling rates were determined from the density, tempera- 
ture, metallicity of the gas and the radiation background. 
While the IGM and ISM are expo sed to UV radiation, 
widely-used metal c ooling rates (e.g., ISuthcrlan d fc Dopital 
1 19931 : iGnat fc Sternb erg 200^) were usually calculated as- 
suming coUisional ionization equilibrium (CIE), which is 
only valid wh en the radiation backgr ound is absent. A re- 
cent study bv lWiersma et al] (|2009al ) has also investigated 
the effects of radiation on cooling in the temperature range 
10^ < T < 10* and also shows that UV largely suppresses 
the cooling rates and shifts the peaks of the cooling curves 
to higher temperatures. 

The radiative cooling was separated into three compo- 
nents: 



A = Ah, He + Ametal + Acon 



(4) 



where Ah, He is net cooling due to primordial species (H, H+, 
He, He"^ and He"^^), Ametai is the rate due to metals, and 
Acomp is the Compton cooling/heating. For the primordial 
gas the ionization, cooling and heating rates are calculated 
directly from the ionization equations with rates matching 
lAbel et all (ll997D . This enables our simulations to capture 
the non-equilibrium cooling of primordial species. Although 
our study focuses mainly on epochs where the IGM is mostly 
in ionization equilibrium, non-equilibrium effects are impor- 
tant at earlier redshifts, especially during reionization, where 
temperature changes of several thousand degrees occur as a 
result of non-equilibrium heating. 

We adapted the relation from Ijimenez et al. I l|2003l ) to 
estimate the helium mass abundance Y as it varies with 
metallicity. For metal mass abundance Z < 0.1, we assume 
Y = Yp + {AY/AZ)Z, where Yp = 0.236 is the helium abun- 
dance for primordial gas, and AY/ AZ = 2.1 is the ratio of 
helium mass to metal mass produced in stars. If the metal 
mass abundance exceeded 0.1 (~ WZq), Y is linearly de- 
creased so that when Z = 1.0 (100% metals), Y — 0. 

The cooling and heating rates due to metals un- 
der the ultraviolet(UV) radiation background were calcu- 
lated using the phot oionization code CLOUDY (version 07.02 
iFerland et al.l l| 19981 )'). Cloudy output assumes that metals 
are in ionization equilibrium, a good approximation when 
extragalactic UV radiation is present. The UV radiation 
field generally varies its function form in space and time. 
We assume a uniform background obtained from extracting 
the CLOUDY (07.02) built-in extragalactic UV background 
which includes radiatio n from both quasars and galaxies 
l|Haardt fc Madaull2005l ) .The radiation field is a function of 
wavelength, and it turns on at redshift z = 8.9 and evolves 
with redshift to z=0. The same UV background was also 
adopted to calculate the cooling due to primordial species. 

Metal cooling depends on redshift because the radia- 
tion background evolves with time. Although our simula- 
tions traced the formation of alpha elements and iron sepa- 
rately, our metal cooling assumed relative solar abundances 
for simplicity. We used the CLOUDY default solar composi- 
tion, which contains the first t hirty elements in the per i- 
odic table (Dat a compiled from [Crevesse fc Sauvall l|l998l ) ; 
iHolwegej l|200ll ): lAllende Prieto et al.l (120011 . 1200^ see Ta- 
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ble 9 of the CLOUDY version 07.02 documentation for de- 
tails). The metal cooling rates were calculated by subtract- 
ing the primordial cooling rates from the cooling rates of 
the total 30 elements, including H and He. We found that 
the net metal cooling rate is linear with metallicity, which 
allows us to tabulate the rates for solar abundance only and 
scale with metallicity (equation [5|. The metal cooling rates 
can thus be written as: 



Ametal{T, p, Z, Z) 







(5) 



where T, p, z and Z are the temperature, density, red- 
shift and metallicity of the gas, respectively. 

The metal cooling rates were tabulated in log(T) (in 
unit of K, from 100 K to lO^K with Alog(T)= 0.05), log(n//) 
(where uh is the hydrogen number density, from 10~® cm"'^ 
to 10^ cm"^ with Alog(nH) 0.1), and redshift z (8.9 ^ z 
^ with grid size 0.1). The cooling and heating rates were 
interpolated separately in log space in the simulation. The 
RMS of relative interpolation errors for cooling and heating 
are less than 1 percent. We assumed the gas was optically 
thin to the ionizing background for all wavelengths, which 
is valid for the IGM. We included metal cooling down to 
100 K which can occur in dense environments or in the IGM 
that is cooled by adiabatic expansion, although the latter is 
generally not not enriched. At low temperature we turned 
off the molecular hydrogen cooling, because H2 in the IGM 
is not in ionization equilibrium and CLOUDY gave unrealisti- 
cally high H2 cooling rates by assuming equilibrium (also see 
[Smith ct al. (2008)). Beyond z > 8.9 cooling was calculated 
with effectively no UV background. 

Figure [1] shows cooling rates normalized by n|f as a 
function of temperature at various metallicities. The black 
curves show the total cooling rates with z=3 radiation back- 
ground, while the red curves show the CIE metal cooling 
(no UV) plus the cooling of H and He under the UV back- 
ground, a model often adopted in cosmological simulations. 
With or without UV, the presence of metals can increase 
cooling rates by up to several orders of magnitude. In the 
temperature range 10* < T < 10* K, radiation ionizes the 
plasma and reduces the number of ions that can be col- 
lisionally excited therefore decreases the cooling rates sig- 
nificantly (e.g., about an order of magnitude at 10^ K in 
Figure [1]). The cooling peak also shifts from 10^ K to 10^ K 
due to the change of ionization states of the metal elements. 
From 100 K to 10* K, the UV background increases the 
number density of free electrons (rie) hence enhancing for- 
bidden line cooling from low ionized species such as C I, C 
II, Si I, Si II and O I. Forbidden line cooling is the domi- 
nant cooling process in enriched IGM at low temperatures. 
At low densities {uh ^ 10~'*cm~^), photo-ionization also 
creates some higher-ionized species with strong magnetic 
dipole transitions (e.g. Ne V and Ne VI). These ions can cool 
the plasma more efficiently, producing a local cooling peak 
around lO^K (Figure[l]). However with an increase in density 
these ions disappear and the major coolants are the less ion- 
ized species mentio ned, consistent with ISM cooling models 
(e.g. IWolfire et al] 2003). Overall, the radiation background 
increases the metal cooling rate substantially between 100 K 
and 10* K over CIE. Figure [T] shows unequivocally that sim- 
ply adding the CIE model of metal cooling on the primordial 
cooling rates results in dramatic over or under estimation of 




Figure 1. Cooling rates normalized by (in erg cm^) with 
metallicities [Fe/H] = -2.0, -1.0, and 0.5. The gas has density 
riff = 10~^ cm~^. Black: total cooling rates (include H, He and 
metals) under the z=3 extragalactic radiation background from 
iHaardt fc Madaul | |2005|'I : Red: Sum of cooling rates of H and He 
under the same radiation background, and metal cooling rates, 
calculated without UV background 



the cooling rates, depending on the temperature. The am- 
plitude of the UV background used here is however too high, 
as implied in the Lyman a flux decrement evolution in Sec- 
tion l3.3l Thus the effect of photo-ionization on metal cooling 
rate may be overestimated. 



2.2 Turbulent Metal Diffusion 

As a Lagrangian particle method, SPH does not include 
any implicit diffusion of scalar quantit ies such as metals. 
IWadslev. Veeravalli fc Couchmaiil (2008) demonstrated that 
this has physically incorrect consequences for even simple 
processes such as convection and Rayleigh- Taylor instabili- 
ties. Including a simple model for turbulent mixing (using 
a diffusion coefficient, D = C Av hspu based on the pair- 
wise velocity, Av, at the resolution scale , /isph and C ~ 0.1) 
IWadslev. Veeravalli fc Couchmaiil (|2008l ) were able to match 
Eulerian grid code results (which must mix due to the neces- 
sary advection estimates). In particular, it became possible 
to generate similar non-radiative galaxy cluster entropy pro- 
files with SPH as with hi gh resolutioii grid c odes. This was a 
majo r discrepancy in th e lFrenk at al.l ( 19991 ) cluster compar- 
ison. iGreif et al.1 l|2009t ) implemented a similar scheme and 
applied it to simulating supernova remnants. As discussed 
previously, galactic outflows should be highly turbulent and 
thus mixing is essential for IGM studies. 

Turbulent mixing models have a long history in envi- 
ronmental and engineering fluid mechanics. GASOLINE now 
incorporates a mo re robust mix ing estimator similar to that 
first proposed by ISmaeorinskvl (tl963i ) for the atmospheric 
boundary layer, 

dA. 



V{DVA), 



D 



(6) 



where A is any scalar, we use the trace-free shear tensor 
for Sij and h is the measurement scale (here ~ hspn)- This 
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choice for Sij results in no diffusion for compressive or purely 
rotating flows. In SPH terms the diffusion expression for a 
scalar Ap on particle p is computed as follows, 

Sij\p = — y^mg(^Jj|g - Vj\p)yp,iWpq, 
Pp 

<3 
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Sij \p 
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-(^ijlp + S'jilp) — (5ij— Trace S\p, 
C \Sij\p\ hp, 
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where the sums are over SPIf neighbours, q, 5ij is the Kro- 
necker delta, W is the SPH kernel function, pq is the density, 
Vpq is the vector separation between particles Vi\q is the par- 
ticle velocity component in direction i, Vp is the gradient 
operator for particle p (operating on the SPH kernel func- 
tion) and Vp,i is the ith component of the resultant vec- 
to r. The difference between this mo del and th e ones used 
in IWadsIev. Veeravalli &: Couchmanl (|200& ') and iGreif et aP 
l|2009l ) is that the diffusion coefficient is calculated accord- 
ing to a turbulent mixing model instead of simply using 
velocity differences. The coefficient depends on the velocity 
shear hence it better models the mixing in shearing flows. 
On the other hand, if there is no shearing motion between 
two phases of fluids (such as a clumpy medium embedded in 
hot wind fluid and the two move with same speed), then no 
turbulent diffusion is added. The Smagorinsky model has 
also been used successfully in other field such as weather 
modeling. 

A coefficient value of order 0.05 — 0.1 is expected from 
turbulence theory (depending on the effective measurement 
scale, h). It was found that a conservative choice of C = 0.05 
was sufficient to match the cluster comparison. This diffu- 
sion was applied to thermal energy and metals in all runs 
except "nmd_40_256" . Diffusion of density results naturally 
from the noise in the velocity field which moves particles 
around. 



3 GLOBAL PROPERTIES OF THE 
SIMULATION 

Before embarking on a detailed examination of the metal 
distribution and evolution in the IGM, it is worth estab- 
lishing the basic properties of the simulated volume with 
respect to the star formation history (SFH), the Lyman-Q 
forest and the evolution of total neutral hydrogen. These 
properties directly refiect the effectiveness of the feedback 
processes and allow the approach to be calibrated with re- 
spect to observations. 



3.1 Star Formation History 

Figure [2] shows the evolution of global star formation rate 
(SFR) as a function of redshif t. The ob s ervati onal data with 
error bars were adapted from iHopkind (|2004l ) and scaled to 
the same IMF and cosmology as used in our simulations. 
The figure shows that all simulations with moderate resolu- 
tion (with or without metal cooling and diffusion) produce 
the SFH that is consistent with observations for z > 0.5. 
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Figure 2. Star formation rate (SFR) density (in unit of Mq 
yr~^ Mpc~^) a s a function of redshift. Data were compiled from 
iHopkinsI ||2004) and scaled to the same IMF and cosmology as in 
our simulations. The moderate resolution simulations give SFR 
consistent with observations down to redshift 0.5. The high resolu- 
tion run has significantly higher SFR but still in broad agreement 
with observations. 



The peak of star formation lies arou nd z ~ 2 for all sim- 
ulatio ns. Without superwind models, Springel fc Hernguistl 
( 20031 1 and lOppenheimer fc Davd (j20od l found higher SFR 
than the observations and the SFH is peaked at z --^ 4. 
Our results suggest that the suppression of SFR can also 
be attained using feedback models without explicit winds. 
Until z ~ 0.5, our 256'^ simulations show no overproduction 
of stars while matching other observables. Note that while 
the 512'^ does show high SFR, it was a direct convergence 
test with all parameters fixed. The parameters affecting SF 
efficiency should be adjusted downwards with increasing res- 
olutio n for consistent outcomes at z = (c.f. IStinson et al.l 
I2006D . 

The simulations do form too many stars at z < 0.5. One 
possible reason is that AGN feedback becomes more relevant 
because large galaxies or clusters become more prevalent at 
lower redshift, and the supermassive black holes in these ob- 
jects can generate strong AGN feedback effects. Since AGN 
were not included in our model, the global SFR is higher 
than observations. 

The simulation without metal cooling produced fewer 
stars at z < 5. The difference in SFR between the simu- 
lations with and without metal cooling increases from ~ 5 
% at z =4 to ~ 20 % at z=l, then decreases to 10 % at 
z=0. Metal cooling enhances the coolin g of the IGM and 
therefore its accretion rate onto galaxies. IChoi fc Nagaming 
(j2009l ) also found enhancement of SFR with metal cooling, 
though their SFR was generally enhanced by 20% to 30% 
over all redshifts below z=15, and in some runs by 50% 



at z=l. The discrepancy may result from their adoption 
of ISutherland fc Dopital (| 19931 ) metal cooling rates that are 
significantly higher than what we used, sinc e it does not 
include the radiation background. However, 'Schayc et alj 
(2010) used the cooling model from Wiersma et al.. (2009^) 
which also includes the suppression of metal cooling by UV 
background, ye t still found a larger enh ancement of SFR 
comparable to IChoi fc Nagamind (|2009l l. This is possibly 
due to the higher metallicities in stars and the ISM in their 
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simulations jWiersma et al.ll2009bl ). In particular, they in- 
cluded mass loss by AGB stars which may largely contribute 
to the enrichment of the ISM. 

The simulation without metal diffusion (but with metal 
cooling) also produced slightly fewer stars (about 6 % mean) 
than the reference run at z < 4, because metal diffusion al- 
lows enriched particles to mix their metals, so more particles 
experience metal cooling and turn into stars. On the other 
hand, the metallicities of the original enriched particles are 
reduced because of the diffusion, thus those particles have 
less strong cooling hence are less likely to form stars. Our 
result suggests that the first factor dominates. Since the ef- 
fect of metal diffusion relies heavily on metal cooling, it is 
smaller than in the case where the metal cooling was turned 
off. 

The effect of resolution on the SFH is shown in the 
convergence runs "mcd_45_256" (cyan) and "mcd_45_512" 
(blue). At z = 8, the high resolution run produces 3 times 
larger SFR. Although the result does not converge by z = 2, 
the difference between these two runs does steadily decrease 
with time. At z = 2, it reaches about 50%. The SFH from 
"mcd_45_512" is close to the high end of the observational 
data but still consistent with them. The increased star for- 
m ation is the sa.me ma gnitude seen in the convergence tests 
m IStinson et ahl (|2006l ). where it was determined that higher 
resolution simulations produce more high density gas. High 
density gas forms more stars either because it surpasses the 
density threshold where gas in the low resolution simulations 
does not, or because star formation is a function of gas den- 
sity in the Kennicutt-Schmidt law. Since the difference is 
most prominent in the early epochs of our simulations when 
galaxies are small, it is likely that many dwarf galaxies sur- 
pass the gas density threshold in the high resolution that 
do not in the low resolution simulations. The effect of initial 
condition variation is also worth noting: the "mcd_45_256" 
run has slightly lower (but comparable) resolution than the 
reference run. However, it has a higher SFR at z > 1. 

3.2 Evolution of Qm 

The evolution of the cosmic mass density in neutral hydro- 
gen (S7Hi)with redshift and its relation to stellar mass den- 
sity is one of the key observables for understanding the in- 
teraction between gas an d galaxies. Observation s of damped 
Ly-g systems (DLAs) (Prochaska et al.l l2005l : iRao et al] 
l2006l : [Prochaska fc Wolfe |2009|) and the H I 21-cm emis- 
sions ijZwaan et al.ll2005l : iLah et al.ll2007l ') (shown in Figure 
O suggest that Qm does not evolve substantially through 
cosmic time, even though the stellar mass density keeps 
increasing. This implies that there is a steady supply of 
gas cooling onto galaxies, providing fuel for star formation. 
The decrease of H I at z ~ 2.3 has been linked to violent 
feedback processes, including SN, galactic winds and AGN 
activities (|Wolfe et al.ll2005l : IProchaska fc Wolfell2009l ). but 
it is not clear why the amount of H I increases back at 
z < 2.3. IProchaska fc Wolfd (|2009l ') argued that the data 
at z < 2 was biased high and JIhi should remain constant 
from z = 2.2 to z—Q (the last 10 Gyrs). The non-evolution 
of 57hi indicates that star formation is self-regulated such 
that gas accretion, star formation and feedback processes 
balance each other. Numerica l simu lations without mass 
loss through winds l|Cen et al.l |2003| ) found Qm is several 



times higher. iNagamine et al.l (l2004a ) i nclud ed the super- 
wind model from ISpringel fc Hernguistl l|2003h in their sim- 
ulations so that neutral gas could be ejected from the galax- 
ies, which produced the correct amount of fiui at z > 2, but 
at z < 2 too much material was blown away, resulting in a 
deficit of Qui (dash line in Figure O . 

Figure |3] shows the evolution of Ohi in our simula- 
tions. Because most neutral hydrogen resides in damped 
Ly-a systems, which are c louds with H I column densities 
larger than 2 x 10^° cm ~^ (jWolfe et al.i 2005. and references 
therein), the gas is mostly self-shielded from external ioniz- 
ing photons. Self-shielding was not modeled during the sim- 
ulations, however we used a radiative transfer post-processor 
to recalculate the ioniza tion states of hydrogen as described 
in iPontzen et al.l ll2008i). The ion izing background was the 
standard Haardt fc Madaul l|2005l ) reduced by a factor of 2 
(as discussed in section [3.31 below). However, variations of 
this magnitude have practically no impact on the dense gas 
that dominates flm- 

Figure |3] shows that after we applied the self-shielding 
correction, our feedback model produced similar flm to ob- 
servations from 2: ~ 3.5 down to z — for the runs with 
moderate resolutions. At higher redshift these simulations 
underestimate ilm, possibly because of insufficient resolu- 
tion (discussed below). The result indicates that our feed- 
back models effectively moderate star formation but not too 
strongly, so that gas accretion onto disks is not disrupted. 
Thus it can maintain the steady supply of neutral hydrogen 
to galactic disks. The shape of the Qui relation follows the 
SFH, which reflects the relation between the H I density and 
SFR density, i.e. the Schmidt Law for star formation. The 
decrement a,t z — 2.3 is not reproduced in our simulations, 
indicating that more violent feedback maybe necessary at 
this limited redshift range if this feature is confirmed. We 
will investigate wind models in high resolution single galaxy 
and galaxy groups in future work. 

As expected, Figure [3] shows that metal cooling and 
metal diffusion increases the amount of neutral hydrogen. 
The metal cooling effect is obvious at z < 3 and at z = the 
increase is ~ 17 %. Metal diffusion has a smaller effect. The 
reasons for these effects are similar to the ones discussed in 
SFR analysis in section |3TT] 

With eight times higher resolution, the simulation 
"mcd_45_512" produces a significantly higher f^ni and a fiat- 
ter curve at z > 2.5 which make the result compare better 
with observations. At z = 5, the high resolution run con- 
tains about 80% more neutral hydrogen. At z = 2 the Qui 
curve joins the moderate resolution run "mcd_45_256" but 
has a slightly steeper decreasing slope towards lower z. The 
difference in f^Hi is likely because the high resolution run 
better resolves small halos, thus increasing the amount of 
self-shielded gas. Since small halos dominate at high red- 
shift, the difference is most obvious there. We note that as a 
result of cosmic variance the "mcd_45_256" simulation con- 
tains 15-30% less HI throughout the simulation than the 
fiducial run even though it has a higher SFR. Thus the re- 
lationship between SFR and the gas reservoir Qui is not 
always monotonic as seen in the cases for "nmc_40_256" and 
"nmd_40_256" . Higher SFR consumes more neutral gas thus 
may decrease the amount of Qui- 
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Figure 3. Evolution of total neutral hydrogen density in units of 
critical density today ( Ohi)- Obser vational da t a poin ts: Square: 



iProchaska et al.l ll2005li ; D iamond: iRao et al.l | |2006|) : Asterisk: 
l|Lah et al.ll2007r) : 7Vmno/e: [Zwaan et al.l l|2005l '). 



3.3 Lyman-a Flux Decrement 

The Lyman-a forest traces the neutral hydrogen in the dif- 
fuse IGM. Since the IGM is optically thin, the mean flux 
decrement of the Lyman-a forest, Dly-q, is sensitive to the 
UV radiation and can be used to calibrate the UV back- 
ground in simulations. We generated the Ly — a spectra 
by tracing a lines of sight through every simulation box. 
Each line of sight consists 25 lines segments, which were 
spread through the entire box such that the periodicity of 
the box was used to eliminate any discontinuity in the result- 
ing spectra. Particles that lay within two smoothing lengths 
measured perpendicularly to their closest line segments con- 
tributed to the absorption along the line, and the amount of 
contribution was weighted using the SPH smoothing kernel. 
The particles' peculiar velocity, thermal velocity and Hub- 
ble expansion were all taken into account. The H I, He I 
and He II fractions were calculated during the simulation 
as described in Section [2.11 The Ly-a optical depth at each 
velocity due to Doppler broadening was calculated and then 
convolved with the Lorentz proflle. 

Figure m d epicts the evolut i on oi Dlvci compared to th e 



observations of [ Bernardi et al. ( 2003^. Schave et al.l 
iKirkman erahl l|2007^ " and iFaucher-Giguere et al.l 



2003), 



200S 



When we attenuate the UV spectrum by factor of 2, our re- 
sults agree with most observations up to z ~ 3.5. At higher 
redshift, our results are slightly lower than, but still broadly 
consistent with the observation al data. The biggest differ- 
ence is seen when comparing to lBernardi et all l|2003l ). The 



discrepancy may indicate that more UV attenuation for z 
> 3 in our simulation is necessary. However, observations of 
the Lya forest at high redshift may be subj ect to large un- 
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) and metal 


The 


Schave et al.l 



data are contamination corrected and our result compare 
better with these data. Hence, a factor of 2 seems sufficient 
and we use this value for all subsequent calculations. This 
attenuated UV background is only used in the analysis, but 
not in the entire simulation. Although for the diffuse, highly 
ionized primordial IGM, the effect of photo-ionization on gas 



Figure 4. The evolution of the mean flux decrement of the Ly- 
a forest from our reference simulation ( "mcd_40_256" ) and the 
convergence test ( "mcd_45_256" and "mcd_45_512" ). Color sym- 
bols wi th error bars indicate the observational data. Green di- 
amond: iBernardi et al. I ll2003h . data here are derived from their 
original evolu tion of effective opt ical depth, Teff in their Fig. 4; 
Blue triangle: ISchave et al. ('2003) , data here are derived from the 
original r^ff after removal of pixels contaminate d by metal lines 
(Table 5); Red solid dot: IKirkman et al] | |2007|) ; Cyan asterisk: 
iFaucher-Giguere et al derived from their Teff evolution 

after contamination correction (Table 1). 



dynamics is small (|Croft et al.lll993) . UV does affect metal 
cooling and hence star formation. Thus, by using the orig- 
inal UV background without attenuation, our metal cool- 
ing calculation may have overestimated the effect of photo- 
ionization. However, along with other workers, we are also 
neglecting local ionizing sources which correct the UV up- 
ward to a significantly greater extent within galaxies. These 
details of ISM physics are potentially important but beyond 
current capabilities where we must rely upon relatively crude 
ISM models. 

Focusing on the convergence test, with high resolution 
(dot-dashed line) the Lya decrement evolves very similarly 
with redshift, but is a few percent smaller than the lower res- 
olution (dashed line) run. The decrease in DLya is consistent 
with the increase in SFR and Qm- It can be interpreted as 
when more gas accretes onto galaxies and forms stars for 
high resolution, the H I content in the IGM (probed by the 
Lya forest) decreases. 



4 EVOLUTION OF GAS, METAL FRACTIONS 
AND METALLICITY 

We now examine the enrichment of the IGM through the 
metal distribution in different gas phases. We examine how 
they evolve and the effect of the metal diffusion and cool- 
ing. The gas phases w ere defined using the convention of 
IWiersma et al.l (|2009lf ) . All the gas that has hydrogen num- 
ber densities uh > 0.1 cm~^ is star forming gas (SF gas). 
Of all the non-star forming gas (non-SF gas), gas with 
lO^K < T < 10'' K is the warm-hot intergalactic medium 
(WHIM) and gas with T > K is the intracluster medium 
(IGM). Of the cooler medium (T < 10^ K), that with over- 
density p/pmean > 100 is associatcd with galactic halos. 
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while that with p/pmean < 100 is the diffuse IGM. Note 
that these definitions are somewhat arbitrary and there is 
no well-defined border for each gas phase. In particular, at 
low redshift substantial part of the WHIM comes from the 
shock-heated gas accreting onto gala^cies so it is likely to 
be associated with galactic halos. However, by this defini- 
tion gas is roughly divided according to its location and key 
physical processes. 

4.1 Evolution of Gas and Metal Fractions 

Figure [5] shows how the gas and metal fractions evolve. The 
left panels focus on the mass fractions in gas and stars. Al- 
most all baryonic mass is in the diffuse IGM at z=7. By z=0, 
the mass fractions in stars, the diffuse IGM and the shock 
heated WHIM are about 20%, 40% and 40%, respectively 
for the moderate resolution runs. The stell ar mass frac- 
tion at z =0 is consistent with observations (|Wilkins et al.l 
|2008| . and references therein), although close to the high 
end. A large amount of WHIM forms at low redshift when 
gas is shock heated when falling into halos. The mass frac- 
tion of WHIM gas is consistent with SPH simulations using 
expli c it galactic superwind m odels (|Oppenheimer fc Pavel 
l2006l: IChoi fc Na.gaminel |2009| ). and Eulerian simulations 
i Cen fc Ostrikeij |2006| ). The amount of star forming gas 
(ISM) and galactic halo gas evolves similarly to the SFH 
shown in Figure [2l reflecting a close relation between the 
ISM, halo gas and SF activity. 

The blue and cyan lines in Figure [5] show the conver- 
gence of mass fractions. At high redshift, the 512'^ simulation 
has a larger mass fraction in stars, the ISM, halo gas and 
the WHIM. Resolving smaller halos at large z increases the 
amount of gas in galax;ies (halos and the ISM) and enhances 
SF. Consequently, enhanced stellar feedback increases the 
amount of WHIM. This effect is more significant at higher 
redshift and decreases with time. The WHIM fraction con- 
verges at z < 4. The ISM and halo gas fractions reach similar 
values as in the low resolution run at z —2, although they 
seem to decline slightly faster towards low z. For the stellar 
mass, however, there is about 37 % more in the high reso- 
lution run at z = 2. Although it is still difficult to predict 
the stellar fraction at z = 0, we expect that the difference 
will be lower at z = since there is already no excess in the 
reservoir gas (the ISM) at z —2. In fact, early star formation 
seen in the high resolution simulation may consume the gas 
that would otherwise contribute to SF at low redshift, hence 
the SFR may decline faster towards low z, which may imply 
a rapid convergence in stellar mass fraction. 

The right panels of Figure [5] show how the metal frac- 
tions present in stars and different gas phases evolve. At z 
— 7, gas (SF and Non-SF) contains the majority of the met- 
als while stars contain only a few percent, which is possibly 
because the time scale for gas consumption is longer than 
the age of the Universe, so that metals in the ISM do not 
have enough time to be incorporated into stars. With metal 
diffusion, initially the non-SF gas (the IGM + halo gas) and 
the SF gas (i.e. the ISM) contain comparable amounts of 
metals, with the former having slightly more in most runs. 
With time, the metal fraction in the non-SF gas decreases 
and the metals budget is dominated by the ISM and later 
by stars. As z=0, 80% to 90% of metals reside in stars. The 
IGM and halo has about 10% of the metals and the ISM 



has only a few percent. The result that the IGM contains a 
larger fraction of metals at higher redshifts suggests that the 
enrichment process is very efficient at early epochs, possibly 
because it is easier for wind material to leave the shallow 
potential wells of early objects. For different phases of the 
IGM, metal fractions in both the WHIM and the diffuse 
IGM decrease with time. For example, in the standard run, 
the WHIM metal fraction decreases from ~ 40% at z=7 to 
about 10% at z=0, the metal fraction in diffuse IGM de- 
creases from ~ 10% to ~ 1%. The IGM, on the other hand, 
increases its share of metals as galaxy clusters form at low 
redshift. 

Metals in the IGM are primarily in the WHIM, with 
a smaller fraction in the diffuse IGM. This remains the 
case for the entire c osmic history, which differs from 
iDave fc Qppenheimeil l|2007l '). who found that most met- 
als reside in the cooler diffuse IGM from 1.5 < z < 6. 
This is probably because they use hydrodynamically de- 
coupled, kinetic feedback. Their superwinds were generated 
in cold SF gas and hydrodynamically decoupled when they 
left the ISM, and the wind material was likely to be cool. 
Our feedback model ejects SN energy into the surroundings 
and suppresses the gas cooling and the winds are generated 
by thermal pressure, which m akes them more likely t o be 
hot and in the WHIM phase. IWiersma et all l|2009bl 'l ob- 
tained a similar result as ours for the IGM metals using a 
local energy injection and hydrodynamica lly coupled wind 
model from iDalla Vecchia fc Schavd (|2008l '). However, met- 
als in their simulations mostly reside in the ISM at redshift 
z > 2 while in our results the IGM generally contains more 
metals at high redshift. While more controlled simulations 
are necessary to explain the difference, it seems that subgrid 
feedback models alter metal enrichment significantly. 

The high resolution simulation increases the metal frac- 
tion in stars while decreasing it in the ISM and Non-SF gas. 
Between different phases of the Non-SF gas, halo gas and the 
diffuse IGM increase their shares of metals, while metal frac- 
tions in the WHIM decreases. This is expected since with 
high resolution SN feedback events are more frequent but 
each has smaller impact. As a result, particles that receive 
feedback energy and metals generally have shorter coolin g 
shut-off times in the blastwave model (|Stinson et al.|[2"006l ). 
Thus they can be cooled more efficiently by metal cooling 
to the diffuse IGM, or to the cool halo gas. No clear conver- 
gence is seen at z = 2 when the high resolution simulation 
was terminated. The differences in stars, the ISM, halo gas 
and the WHIM are 17%, 25%, 14% and 41%, respectively 
However the trend of metal fraction evolution as a func- 
tion of redshift remains similar for high and low resolution 
runs. The general conclusions about enrichment efficiency 
and about metal distribution in each phases also remain the 
same. 



4.2 Evolution of Metallicity 

Figure [6] shows how the mass weighted metallicity in stars 
and the gas phases evolves. Stars and SF gas have the high- 
est metallicities throughout the simulation. The metallici- 
ties are about 10~^ Zq at z — 6-7 and steadily increase 
to ~ 0.5 solar at z=0. At z =0, the value is in agreement 
with th e observed stellar metalli cities in galaxy groups, Z* ~ 
0.6 Zq (|Finoguenov et ani2003l 'l and simulation results from 
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Figure 5. The evolution of baryon mass and metal fraction in various gas phases and stars. The gas is divided as star forming gas, diffuse 
IGM, halo gas and WHIM according to the definition described in the text. Black: the standard run "mcd_40_256" . Red: the run with 
metal cooling turned off. Green: the run with no metal diffusion. Cyan: medium resolution run in the convergence test "mcd_45_256" . 
Blue: high resolution run "mcd_45_5f 2" . The meaning of each line type in each panel is defined by the legends in that panel. 



iDave fc Oppenheimeil ()2007l ). Ho wever, it is about 0.5 dex 
lower than the measurement from lGallazzi et al.l (|2008l ) (di- 
amond symbol in left upper panel of Figure [6]), although our 
total stellar met al density (Qz ) is 4-5 xl0~^ at z = 0, con- 
sistent with the iGallazzi et al.l l|2008l ) value. At z = 2, the 
stellar metallic ity is consisten t , but slightly lower than the 
observation by iHallidav et all |(2008D . 

Metallicities of stars and SF gas are similar throughout 
the simulation. The non-SF gas has an overall similar trend 
of evolution as the stars and SF gas, but is less enriched on 
average. The metallicity of the diffuse IGM evolves signifi- 
cantly but is much less enriched comparing to other phases. 
At z=0 it is only 10"^ '^ Zq. The metallicities of the WHIM 
and the ICM evolve slowly with a slight decreasing trend, 
despite their rapid mass increase at low redshift. Both val- 
ues vary between 0.01 Zq and 0.1 Zq. Since WHIM and 
the cool diffuse IGM contains similar mass at z = 0, the fact 
that WHIM metallicity is much higher than the diffuse gas 



implies our IGM enrichment process happens primarily in 
the WHIM. The convergence is better for metallicity than 
metal fractions, and the metallicities of most gas phases are 
nearly converged at z = 2. The largest difference is seen for 
diffuse IGM and the halo gas, where the high resolution run 
has about 43% and 26% higher metallicites, respectively. It 
may be because feedback gas has a shorter cooling shut-off 
time so it is more likely to cool to < 10^ K, or because 
the high resolution run can better resolve star formation in 
smaller halos at higher redshift, and the metals ejected from 
those objects have not been heated above 10^ K. 

We computed the mean metallicities of the halo gas 
and the ISM and compared it with the metallicity evolution 
of the DLAs and sub-DLAs from Pr ochaska e t al. (200^ 
(triangle symbols in the bottom right panel of Figure |6]). 
The DLA and sub-DLA systems are considered tracing the 
gas in galactic disks and halos. Our results compare well 
with observations in this case. We also compare the diffuse 
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IGM m etallicity with observations usi ng C IV llSchave et al] 
(|2003l '). square symbol) and O VI jAguirre et all l|2008h . 
cross symbol) at z ~ 2.5. Our results are smaller than both 
observations even for the high resolution case. Since the re- 
sults do not converge here, this discrepancy could be due to 
a lack of resolution. Alternatively, it may also suggest ad- 
ditional feedback mechanisms which can enrich diffuse IGM 
are probably necessary here. For the intracluster medium, 
observations found highe r metallicities than our results at 
z=0, about 0.2 -0.5 Zq jAguirre fc Schavd |2007| . and ref- 
erences therein). This inconsistency may be because of the 
AGN feedback absent in the simulations. Moreover, the size 
of our simulation also limits the number of galactic clusters 
that can form (below the cosmic average). Since our simu- 
lations were designed for studying metal enrichment in the 
IGM, we will not further address the intracluster medium. 



4.3 The Effects of Metal Cooling and Metal 
Diffusion 

4.3.1 The Effects of Metal Cooling 

The red curves in the left panels of Figure [S] show the evolu- 
tion of mass fraction without metal cooling. At z = 0, metal 
cooling increases the stellar and halo mass fraction each by ~ 
14%, increases the SF gas fraction by 32%, and decreases the 
WHIM fraction by 7%. It reflects that metals enhance the 
coohng significantly at WHIM temperatures (10^ — 10^ K) 
and thus more WHIM cools onto galaxies and forms stars. 
The mass of diffuse IGM is almost unaffected by metal cool- 
ing (~ 1% decrease at z =0) as its metallicity is very low, 
as shown in Figure [S] 

The right panels of Figure[5]show that the metal cooling 
decreases the metal fraction in the WHIM and increases it 
in the halo gas and the SF gas, because it enhances the 
cooling of enriched WHIM onto galactic halos and disks. 
At z=0, the decrease in WHIM gas is ~ 24%, the increase 
in halo gas is ~ 42%. For stars and SF gas the increments 
are 14% and 32%, respectively. For the WHIM and the halo 
gas, the effects of cooling in metal fractions (which traces the 
enriched gas) are larger than in mass fractions (which traces 
the total gas). For star and SF gas, the effect of metal cooling 
is similar in both mass and metal fractions. For the diffuse 
IGM, metal cooling increases its metal fraction at z > 1.5, 
but decreases it at lower redshift. The early increase is likely 
due to the cooling of WHIM to become diffuse IGM, while 
the later decrease is probably because of the diffuse IGM 
itself is enriched enough so that metal cooling can enhance 
its accretion onto halos. 



4.3.2 The Effects of Metal Diffusion 

Metal diffusion produces similar effects as metal cooling on 
mass fraction. Metal diffusion contaminates large amounts 
of otherwise pristine gas with metals, enabling it to cool 
through metal lines from WHIM to SF gas or halo gas, and 
enhances the SFR. As it depends on the metal cooling, the 
effects of metal diffusion on mass fraction are smaller in 
magnitude than having metal cooling off (similar to the ef- 
fects on the SFH). At z = 0, the increase in halo and stellar 
mass fractions are about 7% and 5%, respectively, and the 
decrease in WHIM is ~ 3%. 



The green curves in the right panels of Figure [5] show 
the direct impact of metal diffusion on the metal content of 
various gas phases. Metal diffusion increases the metal con- 
tent of stars, the ISM and halo gas (i.e., gas in the galaxies), 
and decreases the metal content of the WHIM, the IGM and 
the diffuse IGM (i.e., the IGM). This is counter intuitive 
since one expects diffusion to distribute the metals more 
evenly. However, as the enriched gas is ejected out of the 
galaxies, metal diffusion mixes its metals with the ambient 
gas (i.e., the ISM and the galactic halos) along the outflow 
trajectory so that by the time the outflow reaches the IGM 
(WHIM, diffuse gas or the IGM), its metal content has de- 
creased. In other words, diffusion prevent highly enriched 
gas from transporting all its metals to the WHIM. This is 
further discussed in Section [S] Metal diffusion may alleviate 
the problem found in previous SFH simulations where the 
metals are too inhomogeneous (e.g., Aguirre et al. 200^). We 
will make a detailed analysis of the effects of metal diffusion 
in observable ions such as C HI, C IV and O VI in a future 
paper. 

The impact of metal diffusion changes with redshift. It 
is more significant at high redshift for the WHIM, the halo, 
the ISM and the stars. For example, from z = 7 to z = 
0, the increase in star metal fraction varies from 0.4% to 
580%, in SF gas from 7% to ~ 60%, in halo gas from 30% 
to 300%. The decrease in the WHIM metal fraction varies 
from 0.6% to 60%. Since diffusion arises from velocity shear, 
the diffusion effects are high near winds. The concurrent 
decrease in the WHIM metals and the increase in stellar, 
SF or halo metals imply that the winds take metals from 
galaxies directly to the WHIM. The fact that the diffusion 
impact is higher at high redshift suggests that the winds 
between these phases are more effective at early times. 



5 DISTRIBUTION OF GAS AND METALS IN 
DENSITY AND TEMPERATURE AT Z=0 

In this section, we examine the detailed distribution of mass 
and metals in the density-temperature (p-T) phase diagram 
at z=0. In Figure [71 the green contours indicate the mass 
distribution of gas across the phase diagram while the la- 
beled darker contours indicate mean metallicities for that 
gas. The left panel shows the standard run while the right 
panel has diffusion turned off. Looking at the standard run, 
it can be seen that gas with T < W^K and p < 10 pmean 
follow a power law equation of state (EOS). It is heated 
by photoionization and cooled by adiabatic expansion and 
has generally not participated in star formation or feedback, 
following the standa rd expectation for the diffuse IGM (e.g. 
iHui &: Gnedin|[l997t ). At p > 10 pmean, the gas distribution 
splits into high and low temperature branches. The low tem- 
perature branch extends to star-forming gas in galaxies. Due 
to metal cooling, gas with p > W^pmean can reach temper- 
atures below the atomic hydrogen cooling cut-off at ~ 10* 
K. The WHIM is apparent as the less dense gas with tem- 
peratures above 10^ K. There is less gas around ~ 10^'^ K 
compared to the peaks at lO'^K and 10^'^K due to the peak 
in metal cooling rates at 10"" to lO'' K. 

The WHIM arises from both metal enriched winds and 
pristine virial shocks. Wind material can get out to very low 
densities but with realistic diffusion operating the typical 
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Figure 6. The evolution of metallicity in stars and various gas phases described in the text. The same legend as in Figure [5] is used . 
The solar metallicity is defined as Zq = 0.0127. Obs e rvatio nal data: Diamond: stellar metallicity at z = llGallazzi et al.l|2008f) . 
Asterisk: stellar metall icity at z = 2 from lHallidav et al.l | |2008|) . Triangle s: the metallic i ty evo lution of DLA and sub-DLA systems from 
iProchaska et al.l ||2003|'). Cro ss: the IGM metallicity traced by O VI from lAeuirre et al.l feOOSh : Square: the IGM metallicity traced by C 
IV from lSchave et al.l booal ). The observations were scaled to the same solar abundance as the simulations. Zq = 0.0127 



metallicities are not that high and the cooling times can be 
very long. The most enriched gas is at the highest densities: 
star forming gas in galaxies. The metallicity of this gas can 
be super-solar. The highly enriched gas (near solar metallic- 
ities) is all above p/ pmean ~ 10*. The metallicity decreases 
Eis the wind propagates and mixes moving to lower densities. 
This steady progression with density can be seen in the hot 
gas (T > 10^ K). The WHIM is enriched up to 10"^ Z0 
to lO"'^ Z0, while the diffuse IGM following the power-law 
EOS is barely enriched. The metallicity contours in Figure 
[7] show that metal enriched IGM is hotter than this diffuse 
phase. 

The metallicity distribution on the phase diagram with- 
out metal diffusion in the right panel of Figure [7] is simi- 
lar, but the metals are less evenly distributed. In particular, 
more low density material becomes enriched. For example, 
some gas around T~ IQ'^K and pj pmnan ^ 0.01 has metal- 



licities up to 0.1 solar. On the other hand, relatively dense, 
cool halo gas {p/pmean > 100, T< lO^K) is less enriched, as 
shown in the plot that the horizontal strip of lower metal- 
licities extends to higher densities. As discussed in Section 
14.31 this is because without metal diffusion, the metals are 
locked into the original wind material and must travel with 
it. 

The overall trend in metallicity versus density can be 
seen in Figure [B] The effect of metal diffusion can be clearly 
seen. The reference run produces a positive, near linear rela- 
tion between the density and metallicity in logarithmic space 
for p/pmean > 1. The slope is larger for the under-dense gas. 
The simulation without metal diffusion increases the under- 
dense gas metallicity by a factor of 3 while it decreases it 
for halo gas (with overdensities larger than 100). Overall, 
it gives a shaUower slope in the range of log(p/p 

mean J — 

[0, 3]. The result from the simulation without metal cool- 
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Figure 7. Distribution of Gas and Metallicity in the p-T phase diagram of the reference simulation at z=0. The green contours indicate 
the gas mass distribution, and the blaclt, magenta and blue contours with labels indicate the metallicities. From low to high, the metallicity 
labels are [Fe/H] = -3.0, -2.0, -1.5, -1.0, -0.5 and 0.0. Solar metallicity is defined as Zq = 0.0127. The left panel is the standard run with 
metal diffusion and the right panel is without diffusion. 
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Figure 8. The metallicity-density relation at redshift zero for 
all the simulation runs. The effect of metal diffusion is shown 
in the excess at the low over density range and the decrease of 
metallicity in the halo gas. 



three variables. As seen in the leftmost panel, metal cooling 
is efficient at WHIM temperatures and shifts a substantial 
amount of material from the WIflM to lower temperatures 
and the sharp peak of gas at lO'^K associated with primor- 
dial cooling is removed. Without diffusion, super-enriched 
gas is present which is able to cool to nearly 100 K. With- 
out diffusion, enriched gas is also able to travel to extremely 
low densities as seen in the second panel. 

The impact of metal diffusion is clear in the metallicity 
distribution (the rightmost panel in Figure |9]). With metal 
diffusion, there is a large increase in material at low metal- 
licites, Z < 0.1 Zq). Without diffusion the metals concen- 
trate in a spike just above 0.1 Zq associated with the early 
distribution of metals from a star formation event. Turbu- 
lent mixing spreads the metals from the wind material to the 
surroundings so that substantially more gas is contains met- 
als at low l evels. An incr ease in low metallicity gas was also 
achieved in |w iersma et al. (J009b) using a local smoothing 
technique, but a physically motivated diffusion model results 
in the metals being substantially redistributed in space as 
seen in the second panel of the figure. 



ing is also plotted, although the effect is relatively small. 
Observationally, the metallicity-density relation can be in- 
ferred fro m certain metal tra cers in the QSO spectra. For 
example, ISchave et al.l (|2003l ) used pixel statistics of C IV 
and found the carbon abundance follows a log-linear relation 
with density at the higher redshifts (z — 1.8-4.1). These ob- 
servations may strongly constrain the feedback and wind 
generation mechanism, but uncertainties due to the shape 
of the UV background and the ionization states of the metal 
tracers may change the r esults substantially, as d iscuss ed in 
ISchave et all (|2003| ) and lOppenheimer fc Pavel (|2006l ). We 
defer detailed comparison of specific species with observa- 
tions to a following paper. 

Figure |5] shows the metal-mass-weighted probability 
density functions (PDF) over temperature, density and 
metallicity. Bimodal distributions are seen in the PDF of all 



6 CHARACTERIZING THE WIND 
MECHANISM 

The primary difference between the simulations presented 
here and earlier enrichment studies is the mechanism for 
wind generation. Earlier simulations explicitly added veloc- 
ity to gas in regions where stars form. Here feedback only 
added thermal energy and prevented gas cooling according 
to recipe developed using simulations of isolated Milky Way 
like disks and dwarf galaxies. We have shown in previous sec- 
tions that our adiabatic feedback can cause mass loss from 
galaxies to enrich the IGM. However, since winds are pro- 
duced dynamically as a result of feedback instead of via a 
wind recipe, it is helpful to characterize the wind generation 
mechanisms. In this section, we examine mass loss as a func- 
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Figure 9. The probability density function (PDF) of the metal mass over temperature (left panel), density (center panel) and metallicity 
(right panel) at z=0. Solid lines: the reference simulation; Dashed lines: the simulation without metal cooling; Dot-dashed lines: the 
simulation without metal cooling. 



tion of halo mass for our baseline simulations. The purpose 
is to understand how our current adiabatic feedback gener- 
ates winds that enrich the IGM. Our results do allow scope 
for more vigorous feedback. Additional feedback such as ex- 
plicit superwinds and subgrid AGN models may be included 
in our future enrichment studies. The results presented here 
serve as a starting point for further investigations testing 
more complicated feedback scenarios. 

We identified overdense halos using a friends-of-friends 
(FOF) group finder with linking length e = t (inter-particle 
separation ~ 30 kpc) to find material inside regions of 
overdensity ^ <;125. We compared the baryon fraction 
{Mbary/Mtot) of material inside the halos (black solid line) 
as a function of total halo mass at z = with the cosmic 
mean {Qb/^m, dotted line) in the upper panel of Figure [TOl 
The baryon content is low for low mass halos (< 10^" Mq). 
It increases quickly with mass in the intermediate range, 
and stays near the cosmic mean beyond ~ 10"Mq. Halos 
of masses less than ~ 10^^ contain less than the cosmic 
mean. Using a lower resolution (2 x 128'^ particles) run, we 
ruled out that the decrease below ~ 10^^ Mq is a resolution 
effect. This decrease indicates that gas was either prevented 
from accreting onto those halos, or ejected/stripped from 
them. Possible mechanisms that cause the mass loss include: 
1. Background UV radiation heating the gas and prevent- 
ing it from accreting onto dwarf galaxies. 2. Tidal stripping 
during merger events. 3. Winds generated from the stellar 
feedback mechanism. 

We investigated where the "lost" baryons are at the cur- 
rent epoch by pairing each gas particle with its dark matter 
(DM) partner at the same location in the initial conditions. 
We identified the gas particles that remain in the IGM (i.e. 
not belong to any groups) even though their partner dark 
matter particle has accreted to a halo at z=0. When plot- 
ting any properties of this gas as a function of halo mass, we 
mean the mass of halos to which the DM partner of the gas is 
belong. Note that this method only considers the difference 
between the initial condition and the snapshot at z = 0, so 
it does not capture the information in between such as re- 
accretion of gas after ejection. Hence the result here is used 
as a qualitative indication of the fate of the "lost" gas. A 



more precise description would require track gas particles. 
According to its temperature, this "lost" gas was catego- 
rized as WHIM (lO'^K < T < lO^K ) or cooler, diffuse gas 
(T < lO^K). The bottom panel of Fi gure 1101 shows the mass 
fraction of this gas that is in the form of WHIM (black line 
with error bars). For comparison, the gas inside halos was 
also categorized in the same way as warm-hot halo gas and 
cool halo gas, and the fraction of the halo gas that is in the 
warm-hot form is shown in the same panel (red line with 
error bars). Despite of the large error bars, the curves show 
that gas associated with smaller halos tends to be cooler 
than that associated with large halos. This is expected be- 
cause larger halos have larger gravitational potentials and 
higher virial temperatures. 

We investigated the role of stellar feedback and winds 
in the efficient mass loss by tracing the density history 
of the "lost" gas up to z = 5 to identify winds. If its 
highest overdensity, p{z)max/ pmean, was larger than 100, 
and its current density is less than half of the maximum 
{p{z = 0) < 0.5 p{z)max), then the gas is considered ejected 
as wind material. Since the gas density was measured only 
when snapshots were generated every Az — 0.25, it is pos- 
sible that some wind gas was omitted. However, we ver- 
ified this density method in the non-diffusion simulation, 
where enriched winds can be identified by the metallicities 
of the particles. Our density method successfully detects the 
same material identified by metals. The blue and red lines 
in the upper panel of Figure [10] indicates fractions of the 
halo baryons that were expelled as winds, and are currently 
in the WHIM and cool IGM phase, respectively. We also 
plotted the stellar mass fractions in the same figure (dashed 
line). 

Our results show that baryon loss due to winds is most 
efficient at intermediate mass ranges, 10^° to 1O"M0. This 
is also the range where the stellar fraction increases rapidly, 
which refiects the correlation between star formation and 
wind generation. Below 10^" Af©, baryon fractions of the ha- 
los are very low and there are almost no stars and no winds. 
We traced the accretion history of this gas and found it 
has never accreted onto any objects. Though higher reso- 
lution is necessary to investigate these objects, it is likely 
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that heating due to the background UV radiation prevents 
gas accreting onto such halos, since UV heating can bring 
gas temperature to lO'^K — 10^ A", comparable to the virial 
temperature of dwarf galaxies with mass < lO^'^ Msun- Wind 
fractions start to decrease when M > IO^^'^Mq. This may 
be caused by 1. decrease of wind escaping efficiency with in- 
creasing halo gravitational potentials; or 2. decrease of star 
formation efficiency as halo mass increases, as indicated in 
the decline of stellar fraction in Figure If 01 The latter is 
probably due to the increase of cooling time, as the virial 
temperature of the halo increases and the cooling rate de- 
creases beyond fO*^ K (ref. Figure [T| . For halos with mass 
^ fO^'^MQ, virialized gas can have longer cooling time than 
the dynamical timescale and h ence fuel supplies for S F ac- 
tivities can be largely reduced (|Rees fc OstrikeillT977l ). It is 
worth noting that the stellar fr action from our simulation is 
larger than observations (e.g., IGuo et al] [20091 so there is 
an "overcooling" problem, in which gas in simulations tends 
to cool too rapidly and form much more stars comparing 
with observations. It was suggested that feedback from Ac- 
tive Galactic Nuclei (AGN) may alleviate this problem (e.g. 
[McCarthy et al. 2C)09.) especially for galactic groups. How- 
ever, a discussion of AGN feedback is beyond the scope of 
this paper, and we use the stellar fraction curve only to in- 
dicate SF activities in different mass halos. 

The colored lines in the upper panel of Figure [TO] show 
that the wind material can be in WHIM or cooler IGM 
phases. Galaxies of all masses larger than IO'^'^Mq gener- 
ate winds in the form of WHIM, while only galaxies that in 
range of 10^° — IO^^Mq have cool wind gas. We tracked the 
temperature history for this cool wind gas and found that 
about 60 % of it had T > 10 ' K when it was generated and 
has subsequently cooled. 

Figure [TT] shows the metallicity distribution for the 
baryons within halos, WHIM and diffuse gas as functions 
of halo mass. Note that here the WHIM and the diffuse gas 
shown are those that have their paired DM particles inside 
halos. It is the gas that cannot accrete or has escaped from 
halos. This is a subset of the total WHIM and diffuse gas 
described in Section 3] The metallicity for all gas is low for 
halos with mass < 10"Mq. For the diffuse IGM and the 
WHIM, it is likely because this gas has never accreted onto 
galaxies and experienced star formation. For the gas in ha- 
los, the metallicity is still low due to low star formation. In 
the intermediate mass range, the metallicity of the gas in 
halos increases significantly, because the baryon content of 
halos increases and the galaxies undergo multiple star for- 
mation events and accumulate metals. Consequently, winds 
are enhanced, increasing the metallicities in the diffuse and 
WHIM gas outside those halos. The metallicity of baryons 
in halos saturates above ~ 10^^ Mq, refiecting the decline of 
SF activity with increasing halo mass. The metallicities of 
the escaped material (in both the WHIM and cooler IGM 
phase) decreases, possibly due to the inefficiency of winds 
leaving more massive objects and the decline of star forma- 
tion. Note that the amount of diffuse gas is very small in 
massive objects beyond IO^^Mq so the metallicities have 
large uncertainties. 

The dashed lines in Figure [TT] show the metallicities of 
different gas phases for the non-diffusion run. The effects of 
metal diffusion are seen in two aspects. First, metal diffusion 
contaminates the gas that is prevented from accreting due to 



0.20 




Figure 10. Upper panel: Distribution of baryon and stellar 
mass fractions of halos and wind fraction escaped from halos, as 
a function of total lialo mass in the reference run (z=0). Black 
solid line: baryon fraction (include gas and stars) in halos. Dotted 
line: the cosmic mean value. Dashed line: stellar mass fraction. To 
indicate how halos of different mass lose gas and the role of stellar 
feedback and winds, we identify the lost gas by paring each gas 
particle with its dark matter neighbour in the initial condition, 
and finding the ones who have their DM partners within halos but 
themselves in the IGM. When plotting this gas as a function of 
halo mass, we mean the mass of halos that its pairing dark matter 
particle is belong to. The gas found in this method is further 
split to the WHIM phase (10^ K < T < lO'^K) and the cooler 
IGM phase (T < 10^ K). Among this gas, winds were identified by 
tracing the density history of the gas using p{z)max / pmean > 100 
and p(z = 0) < 0.5 p(z) max . Blue solid line: mass fraction of 
the halo baryons that escaped as wind and are currently in the 
WHIM. Red solid line: mass fraction of the halo baryons that 
escaped as wind and are currently in the diffuse IGM. Lower 
panel: For the gas that is currently within or outside of a certain 
halo, the mass fraction of this gas that is in the warm-hot phase, 
as a function of the halo mass. Black curve with error bars: mass 
fraction of the gas outside of halos (identified by the partner DM 
particles) that is in the warm- hot phase (i.e. WHIM). Red curve 
with error bars: mass fraction of halo gas that is in the warm-hot 
phase. 
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Figure 11. The distribution of metallicities of baryons within 
halos and different phases of gas that cannot accrete onto or es- 
caped from halos (described in the caption of Figure IIOII as a 
function of total group mass at z = 0. The solar metallicity used 
in this plot is Zq = 0.0127. 



UV or stripped by tidal forces, as indicated in the increases 
of the metaUicity of the gas associated with low mass halos 
when metal diffusion is on. Second, metal diffusion allows 
winds to lose metals along the wind trajectory, therefore 
significantly decreasing the metallicity in winds, while in- 
creasing it in halo gas. This is reflected in the metallicity 
decrease of the gas outside halos (both the WHIM and the 
diffuse gas) at mass larger than and the increase 

in the metallicity of the baryons within halos. Without dif- 
fusion almost all metals in the unbound gas are locked in 
wind material (98% of all metals in the WHIM and 99% 
of all metals in the diffuse IGM). In the full model, includ- 
ing metal diffusion, the numbers decrease to 58% and 80%, 
respectively. 



7 SUMMARY AND CONCLUSIONS 

We investigated the enrichment of the intergalactic medium 
with SPH cosmological simulations using an adiabatic stel- 
lar feedback model. The simulations incorporated a self- 
consistent metal cooling model with an ultraviolet (UV) 
ionizing background along with metal diffusion that models 
the turbulent mixing in the IGM and the ISM. It was found 
that the UV background significantly alters the metal cool- 
ing rates at all temperatures from 100 K to 10^ K. Above 10^ 
K it decreases the cooling rate and shifts the cooling peak 
to higher temperature, while below 10* K the UV increases 
the metal cooling rates due to the increase of free electrons. 

The simulations produced an SFH broadly consistent 
with observations to redshift z ~ 0.5, and a steady cosmic 
total neutral hydrogen fraction (SIhi) that compare rela- 
tively well with observations, although possible discrepan- 
cies in the observed fini at z ~ 2-3 allow for more vigorous 
mass-loss. This demonstrates that adiabatic feedback can 
moderate SF while maintaining a regular supply of H I. The 
evolution of the mean flux decrement in the Ly-a forest in 
our simulations is consistent with observations to z ~ 3-4, if 



the m agnitude of the UV backg round is lowered with respect 
to the lHaardt fc Madaul l|2005l ') rates in cloudy. 

As the universe evolves, there is a rapid increase in 
the amount of warm-hot intergalactic medium (WHIM) 
and a decrease in the cooler diffuse IGM. At z=0, about 
40% of the mass is in WHIM, co nsistent previou s sim- 
ulations with differe nt methods (|Cen fc Ostrikeij l2006l : 

lOppenheimer fc Dave|[200^ '). The metal content of the Uni- 
verse evolves from the largest fraction being in the IGM to 
the majority residing in star forming gas, to ultimately being 
locked in stars at the present day. These trends reflect more 
effective wind escape at high redshift. IGM metals primar - 
ily r eside in the WHIM, urr l ikelOppenheimer fc Pavel (|2006l ') 
and iDave fc Qppenheimeil (|2007l ). whose metals largely re- 
side in cool gas, as the wind gas starts cool with kinetic 
feedback and less likely to be heated up when the wind ma- 
terial is decoupl ed from the ISM. Our result is however in 
agreement with |w lersma et al.l (|2009bl ') who also used ki- 
netic feedback but with non-decoupled wind models. The 
mean metallicities of stars, star forming gas, galactic halo 
gas and the cold diffuse IGM all increase with time, but 
those of the WHIM and the IGM remain mostly constant 
with a slight decreasing trend. The metallicity of the WHIM 
stays between 0.01 to 0.1 Zq. The metallicity of the IGM 
is similar to the WHIM, which is smaller t han the observed 
value, 0.2-0.5 Zq (jAguirre fc Schave''2007'). possibly due to 
absent AGN feedback. The metallicity evolution of the gas 
in galaxies c ompares well with DL A and sub-DLA obser- 
vations from iProchaska et al.l (|2003l ) and the metallicity of 
the diffuse IGM at z 2.5 is less than the ob servations 
from lSchave et all (|2003l ') and lAguirre et al.l l|2008ll , suggest- 
ing higher resolution, or additional mechanisms for enriching 
the diffuse IGM are probably necessary. 

We characterized our galactic mass-loss and wind gen- 
eration. For the current adiabatic feedback model, winds are 
most efficient for galaxies in the intermediate mass range of 
IO^^Mq to about 10"Mq. Below 10^°Mq gas is likely to be 
prevented from accreting due to UV heating and remains as 
low- metallicity gas. Above about 10^^ Mq fraction of wind 
gas decreases, possibly because the wind escape efficiency 
decreases with increasing halo potentials, or because the de- 
cline of star formation activities especially for massive halos. 
Most winds were hot when generated, but the ones expelled 
from intermediate mass range galaxies, having temperatures 
~ 10^ — lO^'K, could cool through metal lines and become 
diffuse IGM rather than WHIM. 

We investigated the effect of metal cooling and diffu- 
sion on the SFH, the evolution of ilni and the evolution of 
mass and metals in different phases. For metal diffusion, we 
further studied its effect in the density-temperature phase 
diagram at z =0 and the distribution of baryons as function 
of the galaxy mass. Metals significantly enhance the cool- 
ing of the WHIM, allowing the gas to cool and join galactic 
disks. Metals also enable cooling below 10*K. Metal cooling 
decreases the mass and metal fractions of the WHIM while 
increasing the metals in stars, halos and SF gas and increas- 
ing the SFR by 20% and ^hi by 17% at z =0. With realistic 
diffusion included, metals mix between winds and surround- 
ing gas before they leave the galaxies, decreasing the metal 
content in the WHIM and diffuse IGM but increasing it 
in the galactic halo and star forming gas. It prevents en- 
riched, hot winds from creating highly-enriched low density 
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regions, and makes the density-metallicity relation smoother 
so it follows a nearly log-linear relation in the density range 

log{p/pmean) =0-4. 

We performed two additional simulations with one of 
which having 8 times more particles, to investigate the con- 
vergence of our results and the resolution dependence of our 
model. Full convergence is not seen in most of our results, 
however the results from the high resolution run in SFH, ilui 
evolution, Lya decrement evolution are still broadly consis- 
tent with observational data, and the differences caused by 
resolution decrease with redshift. Also, the basic results of 
metal fraction evolution trend and enrichment efficiency in 
different gas phases remain unchanged. The impact of res- 
olution can be seen in three aspects. Firstly, with high res- 
olution the amount of neutral gas (H I) and SFR increase, 
hence the stellar mass fraction and the WHIM mass frac- 
tion related to wind activities. The difference is largest at 
high redshift. Secondly, the metal fraction locked up in stars 
increases with resolution while decreasing that in the ISM 
and IGM in general, suggesting a decline in wind efficiency. 
Thirdly, with high resolution the diffuse IGM and halos have 
a higher metal fraction at the cost of the WHIM metals. 

In future work, we will make a more detailed compar- 
ison with metal absorption line observations. We will also 
examine the detailed behaviour of various wind models and 
the effect on the immediate environment of galaxies. 
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